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Abstract 

We perform molecular dynamics simulations to investigate the relationship between structural 
order and water-like dynamic and thermodynamic anomalies in spherically-symmetric potentials 
having either one or two characteristic length scales. Structural order is characterized by transla- 
tional and orientational order parameters. We find that (i) dynamic and thermodynamic anomalies 
exist for both one-scale and two-scale ramp potentials, and (ii) water-like structural order anomalies 
exist only for the two-scale ramp potential. Our findings suggest that the water-like relationship 
between structural order and anomalies is related to the presence of two different length scales in 
the potential. 
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Most liquids become denser upon cooling and more viscous upon compression. However, 
water and many other liquids with local tetrahedral order Q (e.g., silica, silicon, carbon, and 
phosphorous) show a decrease in density upon cooling (density anomaly) and an increase of 
diffusivity upon pressurizing (diffusion anomaly). These liquids share many other thermo- 
dynamical anomalous properties than those mentioned above. For instance, experiments in 
phosphorous indicate the presence of a liquid-liquid phase transition |2| and similar results 
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are obtained from computer simulations in silica [3J, silicon |4J, and water pj. A possible ex- 
planation of these anomalies is the tendency of these liquids to form bonds resulting in local 
open structures not present in simple liquids. Therefore, much effort has been expended to 
understand the relationship between the structure and the dynamical and thermodynamic 
anomalies observed in tetrahedral liquids. 

Simple spherically-symmetric potentials with a 'core-softened' or a repulsive interaction 
at short distances are able to generate water-like density and diffusion anomalies 
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and can even show a liquid-liquid transition [12]. These findings in simple 



liquids imply that strong orientational interactions (as observed in water or silica) are not 
a necessary condition for a liquid to have thermodynamic and dynamical anomalies. It is 
still not clear whether such strong orientational interactions are necessary for a liquid to 
have structural anomalies. Here we address the question of whether spherically-symmetric 
potentials are also able to reproduce the structural anomalies found in systems with local 
tetrahedral order. 

Errington and Debenedetti Q| (ED) have studied the microscopic structural order in liq- 
uid water by using simple geometrical metrics or order parameters. The structural order was 
characterized using two different metrics: a translational order parameter t [3], quantifying 
the tendency of particle pairs to adopt preferential separations, and a bond-orientational 
order parameter q quantifying the extent to which a molecule and its four nearest 

neighbors arrange in a tetrahedral local structure, as is the case in hexagonal ice. A useful 

way of investigating structural order in liquids is to map state points onto the t — q plane. 
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Such a representation was introduced by Torquato and coworkers [17], who applied it to 
sphere packings and referred to it as an order map. ED used the order map to investigate 
structural order in water jl^j]. Because of the distinctive features discovered in that study, 
in what follows we refer to water-like order maps as the ED order map. ED found that the 
state points accessible to the liquid state in the order map fall into a two-dimensional area, 
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meaning that in general t and q are independent. However, for those state points where the 
dynamical and thermodynamical anomalies occur, the t — q parameters fall on a line in the 
ED order map, meaning that they are strictly correlated. This is a clear evidence of the 
relationship between structure and water anomalies. Shell et al. tested the ED order 
map in silica and found it qualitatively similar to water. However, t — q points corresponding 



3ut fall on a stripe region. 



to anomalies do not fall on a single line as in water 

The ED order map was studied for hard spheres 
In these systems, where no dynamical or thermodynamic anomalies are observed, the liquid 
state points always fall on a single line on the ED order map. In other words, in contrast 
to the case of silica or water, the q — t parameters are always strongly correlated. In 
light of these findings, it is natural to inquire about the ED order map of systems that are 
spherically symmetric and, in addition, exhibit water-like anomalies in their thermodynamic 
and transport properties. 

Here we perform discrete molecular dynamics simulations to study the model introduced 
by Jagla |g] (see Fig. [TJ) . We will identify the model with A = <ti/<j = 1.76 as the two-scale 
ramp potential (2SRP) model, and the model with <7o = as the one-scale ramp potential 
(1SRP) model. We use NVT ensemble for a system composed by N = 1728 (2SRP) or 850 
(1SRP) particles with periodic boundary conditions and control the temperature with the 
Berendsen thermostat. The details of the simulation are given in Ref. [l3j. However, we 
note that we use different units than in Ref. jl^: lengths are measured in units of a% and 
energy is measured in units of U%. We also check that the results do not depend on the 



number of particles and the value of o\ and U±a> 



We use the translational order parameter 
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ter renormalization. 
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t= / \g(s)-l\ds. (1) 



o 



1 /Q 

Here s=rp„ is the radial distance scaled by the mean interparticle distance, g(s) is the pair 
correlation function, and s c a numerical cutoff that can be set to a suitable value (we choose 
s c so that it corresponds to one-half the simulation box size). For a completely uncorrelated 
system, g(s) = 1, and thus t = 0. For systems with long-range order, the modulations in 
g(s) persist over large distances, causing t to grow. 

An orientational order parameter introduced by Steinhardt et al. tecfl and used in 
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171 Il9t |21[ is modified to characterize the average local order of the system |22 1 . 



For each particle, we define 12 vectors (or bonds) connecting the central particle with each 
of its 12 nearest neighbors. Each bond is characterized by two angles (9, ip) and the corre- 
sponding spherical harmonic Yi m (9, ip) can be computed. The orientational order parameter 
associated with each particle i is 
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Here, Y irn (9, ip) denotes the average of Yg m (9, ip) over the 12 bonds associated with particle 



For i = 6 



15j], Qe has maximum value for most crystals such as fee, hep and bec |20[. 



The values of for each molecule in the system obey a Gaussian distribution, and the 



averaged value of Q 6i over all the particles 



22j, Q 6s , is used to characterize the local order 



of the system. In general, Qq s grows in value as the crystallinity of a system increases. For 
example, in fee lattice, = 0.574, and for uncorrelated system, Q 6s = l/\/l2 ~ 0.289. 
Thus, Qq s provides a measure of orientational order in the system. We note that fee is the 
structure of the stable crystal at low pressure in the 1SRP and 2SRP models 

Figs. |2ta) and|l[b) show the density dependence of t for both the 1SRP and 2SRP at 
different temperatures T. The behavior of t(p) is qualitatively the same in both models. 
At low-T, t(p) shows a clear minimum and maximum and, hence, a range of p where t 
decreases with increasing p. This anomaly implies that the system becomes less ordered 
upon compression. As T increases, the extrema in t(p) disappear. At high-T, for the 2SRP, 
t is a monotonically increasing function of p, as is the case in normal liquids. However, for the 
1SRP, t decreases with p. This behavior at high T is probably a consequence of the absence 
of a hard core in the pair interaction potential. Therefore, upon compression particles can 
pass through each other and, as p increases, the structure of the liquid resembles more the 
structure of the gas. As expected, t(T) and Qq s (T) in both models decreases with T at fixed 
p, meaning that the system becomes more disordered upon heating (Fig.|2J). 

Qq s as a function of density for both models is shown in Figs. Efc) and|2fd). At low 
T, Qq s (p) for the 1SRP is a monotonically increasing function of p, i.e., orientational order 
increases upon compression. Instead, Q§ s {p) for the 2SRP at low T shows a clear maximum, 
indicating a range of density for which Qq s (p) decreases with p. We also note that at 
high pressure/density, the two models have different crystal structures, hep for 1SRP and 
rhombohedral for 2SRP. Thus the different Q& s {p) in liquid state may be related to the 
difference in crystal formation. 



Fig. El shows the isotherms from Fig. El in the ED order map representation. For both 
the 1SRP and 2SRP, the state points fall on a two-dimensional region, i.e., t and Qq s can 



be changed independently. As is the case of silica and water [lj, [18[, we also find an 
inaccessible region in the ED order map where no liquid state points can be found. The 
structural anomalies correspond to the section of the isotherms in Fig. |2] where both Q 6s 
and t decrease upon compression. Only for the 2SRP we find such anomalies, as has been 
observed in both silica and water. The state points corresponding to the structural anomalies 
do not fall, strictly speaking, on a single line delimiting the inaccessible region in the ED 
order map (Fig. Efb)). Instead, these state points form a narrow stripe region, resembling 
water ED order map much more closely than ED order map found for silica {18]. 

Next, we discuss the regions in the phase diagram where the structural, dynamical, 
and thermodynamic anomalies occur. In water, the region in the T — p or P — T plane 
corresponding to both structural anomalies contains the region corresponding to the diffusion 
anomaly which in its turn contains the region corresponding to the density anomaly. In silica, 
the diffusion anomaly region contains the structural anomaly region, which also contains 
the density anomaly region. Fig. H] shows the temperature of maximum density (TMD) line 
and the diffusivity maxima/minima (DM) line. The TMD defines the boundary where the 
density anomaly occurs while the DM defines the boundary where the diffusion anomaly 
occurs. By definition, the structural anomaly region in the T — p plane is delimited by the 
location in the T — p plane of the Qq s maximum and the t minimum. Fig. Efb) shows that 
the relation among the regions of various anomalies for the 2SRP is the same as in the case 
of water i.e., the structural anomaly region contains the diffusion anomaly region which 
also contains the density anomaly region. 

For the 1SRP, there is no clear Qq s maximum (i.e., Qq s shows no anomaly), so we are not 
able to identify a structural anomaly region. However, in this case the structural anomaly 
can be identified by those state points along an isotherm where t decreases upon compression. 
Therefore, the maximum and minimum values of t in Fig. Efa) define the structural anomaly 
region. Fig. a) shows that with this definition the structural anomaly region contains 
the diffusion anomaly region which also contains the density anomaly region. We show in 
Fig. Efb) that for the 2SRP, defining the structural anomaly boundaries using the extrema 
of t does not alter the relationship between the regions of various anomalies. Furthermore, 
comparing Figs. Efa) andEJb), we observe that the effect of reducing the hard core distance 
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o"o is to open the structural anomaly region (curves C and A). 

In summary, we find that the 2SRP shows not only density and diffusion anomalies but 
also the same structural anomalies found in tetrahedral liquids such as silica and water. Fur- 
thermore, the 2SRP also shows the same relation among structural, dynamical and density 
anomalies. Our finding suggests that the water-like relationship between structural order 
and anomalies may be due to the presence of two different repulsive length scales. When 
eliminating the hard core interaction, we find no water-like relation between structure and 
dynamics. This suggests that the ratio between the two length scales in the 2SRP, 01/00, is 
the relevant variable in the interaction potential. 

This work is an outgrowth of fruitful interactions with S. Chatterjee and P. G. 
Debenedetti. We thank P. G. Debenedetti for reading the manuscript and many critical 
suggestions. We also thank NSF grants CHE 0096892 and CHE 0404699 for support and 
Yeshiva University for location of CPU time. 



[1] C.A. Angell et al, Phys. Chem. Chem. Phys. 2, 7549 (2000). 

[2] Y. Katayama et al, Nature 403, 170 (2000). 

[3] I. Saika-Voivod et al, Phys. Rev. E 63, 011202 (2001). 

[4] S. Sastry and C. A. Angell, Nature Materials 2, 739 (2003). 

[5] P. H. Poole et al, Nature 360, 324 (1992). 

[6] P. G. Debenedetti et al, J. Phys. Chem. 95, 4540 (1991). 

[7] F. H. Stillinger and D. K. Stillinger, Physica A 244, 358 (1997). 

[8] E. A. Jagla, J. Chem. Phys. Ill, 8980 (1999); Phys. Rev. E 63, 061509 (2001). 

[9] M. R. Sadr-Lahijany et al, Phys. Rev. Lett. 81, 4895 (1998); Phys. Rev. E 60, 6714 (1999). 

[10] A. Scala et al., J. Stat. Phys. 100, 97 (2000); Phys. Rev. E 63, 041202 (2001). 

[11] P. C. Hemmer and G. Stell, Phys. Rev. Lett. 24, 1284 (1970). 

[12] G. Franzese et al, Nature 409, 692 (2001); Phys. Rev. E 66, 051206 (2002). 

[13] P. Kumar et al, Phys. Rev. E (in press) cond-mat/0411274 

[14] J. R. Errington and P. G. Debenedetti, Nature 409, 318 (2001). 

[15] J. R. Errington et al, J. Chem. Phys. 118, 2256 (2003). 

[16] P. L. Chau and A. J. Hardwick, Mol. Phys. 93, 511 (1998). 



6 



[17] S. Torquato et al, Phys. Rev. Lett. 84, 2064 (2000). 
[18] M. S. Shell et al, Phys. Rev. E. 66, 011202 (2002). 
[19] T. M. Truskett et al, Phys. Rev. E. 62, 993 (2000). 
[20] P. J. Steinhardt et al, Phys. Rev. B. 28, 784 (1983). 
[21] A. Huerta et al^L Chem. Phys. 120, 1506 (2004). 



[22] In Refs. 



15 



Yf. m {9,ip) denotes the average over all the bonds in the system, and the 



orientational order is defined for the whole system, and not just for single molecules. However, 

1^ | the orientational order is defined for each molecule and its average 



for water and silica 
is used. 
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FIG. 1: The ramp potential introduced by Jagla 8j. do corresponds to the hard-core distance, 
o"i characterizes a softer repulsion range that can be overcome at high pressure. The 1SRP has 
o"o = 0, while the 2SRP has two length scales with cti/ctq = 1.76. 
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FIG. 2: The upper panels show the density-dependence of the translational order parameter t. 
The solid lines are polynomial fits to the data, introduced as a guide to the eye. (a) One-scale 
ramp potential (ao = 0). From top to bottom isotherms correspond to T =0.04, 0.05, 0.06, 0.07, 
0.08, 0.09, 0.10, 0.11, 0.12, and 0.13. (b) Two-scale ramp potential with cri/cr = 1.76. From top 
to bottom isotherms correspond to T =0.027, 0.036, 0.045, 0.063, 0.082, 0.109, 0.145, 0.172, 0.2, 
0.236, and 0.290. The lower panels (c)-(d) show the density-dependence of the orientational order 
parameter Q§ s for same sets of isotherms. 
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FIG. 3: (a) The ED order map for the one-scale ramp liquid. Isotherms correspond to T =0.04, 
0.05, 0.06, 0.07, 0.08, 0.09, 0.10, 0.11, 0.12, and 0.13 (from top to bottom), and arrow indicates the 
direction of increasing density. For T < 0.08, t has a maximum at a low density and a minimum 
at a high density but Qq s shows no maximum, (b) The ED order map for the two-scale ramp 
liquid. Isotherms correspond to T =0.027, 0.036, 0.045, 0.063, 0.082, 0.109, 0.145, 0.172, 0.2, 
0.236, and 0.290 (from top to bottom). Q§ s and t have maxima at low densities, and t has minima 
at high densities. We can identify a structurally anomalous region bounded by loci of maximum 
orientational order (at low densities) and minimum translational order (at high densities) in which 
both t and Qq s decrease upon compression. All state points within the structurally anomalous 
region fall in a narrow stripe region adjacent to the inaccessible region where no liquid state point 
can be found. 
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FIG. 4: Temperature of maximum density (TMD) and diffusivity minima and maxima (DM) lines 
for the two potentials in the P — T plane. The region of diffusion anomaly is defined by the loci of 
DM inside which the diffusivity increases with density. The thermodynamically anomalous region 
is defined by TMD line, inside which the density increases when the system is heated at constant 
pressure. 
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FIG. 5: Relationship between the structural order and the density and diffusion anomalies, (a) 
For the one-scale ramp liquid, the open region bounded by the loci of t minima (curve A) and t 
maxima (curve C) [see Fig. Of a)] defines the structurally anomalous region. This region contains 
the diffusion anomaly (delimited by the DM) and the density anomaly (delimited by the TMD). 
(b) For the two-scale ramp liquid, the structurally anomalous region can be defined by the region 
between the loci of t minima (curve A) and by the loci of either t maxima (curve C) or Qq s maxima 
(curve B) [see Fig. Efb)]. The structurally anomalous region contains the diffusion anomaly region 
which also contains the density anomaly region. Only the two-scale ramp potential shows the same 
relation between the structural, diffusion, and density anomaly regions as observed in water. 
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